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We discuss a new approach for the numerical evaluation of loop integrals. The fully 
numerical calculations of an infrared one-loop vertex and a box diagram are demon- 
strated. To perform these calculations, we apply an extrapolation method based on 
the e-algorithm. In our approach, the super high precision control in the numerical 
manipulation is essential to handle the infrared singularity. We adopt a multi-precision 
library named HMLib for the precision control in the calculations. 

1 Introduction 

Measurements of fundamental parameters with high precision will be one of the important 
issues at the ILC experiment. We suppose to determine the masses and couplings for the 
standard model (SM) and some of the beyond the standard models like the minimal super- 
symmetric modcl(MSSM) within a few percent precision pQ. These precision measurements 
require the knowledge of high precision theoretical predictions, especially the computation of 
higher loop corrections. Since, the vast number of Feynman diagrams appear in the loop cal- 
culations, performing such computation is absolutely beyond the human power if it should 
be done by hand. The procedure of a perturbation calculation is well established, thus 
computers must be able to take the place of human hand. In this purpose several groups 
have developed computer programs which generate Feynman diagrams and calculate cross 
sections automatically, like GRACE [2], FeynArts-FormCalc [3] and CompHEP @J. 

The GRACE system, developed by Minamitateya group, is one of the systems to calculate 
Feynman amplitudes including loop diagrams. Our final goal of the GRACE system is to 
construct the fully automatic computation system of multi-loop integrals. In this stage, 
it is successfully working for one-loop calculations in both the SM [5] and the MSSM [6]. 
However, the fully automatic way to compute multi-loop integrals is not established, thus 
the system is still semi-automatic. Especially, since infrared singularities require the analytic 
manipulations, it makes difficult to handle the singularity in the automatic way. To avoid 
such a difficulty, we really need to treat loop integrals in a fully numerical way. 

In this paper, we discuss the fully numerical computation of one-loop integrals with ex- 
isting infrared singularities. First, we demonstrate the calculations of an one-loop vertex 
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and a box diagram in QED case. We put a fictitious photon mass as a regulator of the in- 
frared singularities. For the calculations, we adopt the very brute force way as the numerical 
extrapolation method. In our approach, the e-algorithm is efficient. We also introduce a 
precision control technique to handle the numerical calculation with extremely small photon 
masses. 

In QCD, one can not avoid to use the dimensional regularization to treat the infrared 
singularity. The sector decomposition technique [7] is useful to pick up the singularities but 
in general, the numerical integration of each term is not trivial because the denominators 
of the integrands can be zero in the integral region. Therefore, the extrapolation method 
should be effective in this case, too. 

The layout of this paper is as follows. In section 2, we show the basic idea and formula 
of the loop integrals. Numerical results of the infrared one-loop vertex and box diagram are 
shown in section 3. The case for the dimensional regularization is discussed in section 4. In 
section 5 wc will summarize this paper. 



2 Basic idea 



The basic idea of our numerical method is the combination of an efficient multi-dimensional 
integration routine and an extrapolation method. Here, we consider a scalar one-loop n-point 
integral given by 



dH 



(2tt) 4 z (P-m\ + ie){{l+ Pl ) 



ie)' 



(i) 



where I is the loop momentum, pj is the momentum of the jth 
external particle and rtij is the mass carried by the jth internal 
line. Putting e > in the denominator of Feynman integrals to 
prevent the integral from diverging, we can get the numerical 
results of /(e) for a given e. Calculating a sequence of I(ei) 
with various value of ei = eo x [const.)~ l (I = 0, 1, 2, • ■ ■ ) and 
extrapolating into the limit of e — *■ 0, we can get the final result 
of the integration which appears in a physical amplitude. 

In this approach, we use DQAGE routine in the multi- 
dimensional integration. It is included in the package 
QUADPACK [8j and it is a globally adaptive integration routine. 
We also apply the e-algorithm introduced by P.WynnJS]. The 
algorithm accelerates the conversion of the sequences. Our 
approach is very efficient for massive one-loop and two-loop 
integrals including non-scalar type integrals [TU1 HH HH H3] • 
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diagram 



One-loop vertex 



A [GeV] 


Average lost-bits 


Maximum lost-bits 


10 -.u 


88 


92 




98 


102 


1Q -22 


108 


112 



Table 1: Lost-bits information from HMLib with the quadruple precision arithmetic. A is a 
fictitious photon mass. 
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In the above massive cases, the results are consistent with analytic ones very well, even 
in the double precision arithmetic or quadruple precision arithmetic. On the other hand, for 
massless cases, we need much more careful treatments to handle the infrared singularities. 
For the infrared divergent diagram, it becomes harder to get a result with an enough accuracy 
even in the quadruple precision arithmetic. Thus, the precision control of the numerical 
integration becomes another important elements in our approach |14j . 

For the precision control, we use HMlibrTS] as the multi- 
precision library. The advantage of this library is that it gives 
an information of the lost-bits during the calculation, thus we 
can guarantee the precision of the results. In Tabled! we show 
an example of the lost-bits information supplied by HMLib. 
Here, we demonstrated the calculation of one-loop vertex di- 
agram (Fig[T]). The diagram includes the infrared singularity 
due to the photon exchange. We put a fictitious photon mass A 
as a regulator of the infrared singularity. In P-precision presen- 
tation implemented in HMlib, a sign bit is 1 bit and an exponent 
bit is 15 bits and a mantissa is (32 x P — 16) bits with based 
on IEEE754 FP. When P = 4 (quadruple precision), the man- 
tissa becomes 112 bits. In the Table when A is 10~ 23 GcV, the 
maximum number of lost-bits is the same as the number of bits 
of the mantissa. This means quadruple-precision is not enough 
and we have to perform the calculation within octuple or higher precision. 
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Figure 2: One-loop box di- 
agram in QED 



3 Numerical results 



For the one-loop vertex diagram shown in Fig. [TJ the loop integral we consider in this paper 
is 



where 



D 



1 j-l-X 

dx / dy- 
Jo u 

\2 2 



xys + (x + y) m + (1 — x — y)X . 



(3) 



Here s denotes squared central mass energy and m and A are a mass of external particles and 
a fictitious photon mass respectively. Replacing s by s + it in |3]) and applying extrapolation 
method, we obtain the numerical results which are shown in Table [5] [TJj. Results are 
compared to analytic results evaluated by the formula in [16j . The combination of the 
extrapolation method and multi precision control with HMLib works very well and we get 
stable results even though the photon mass becomes much smaller as 10~ 160 GeV. 

We also demonstrate the calculation of the one-loop box diagram shown in FigO Here, 
we consider the following integral, 



dx 



dy 



(4) 



where 



D = —xys — tz(l — x — y — z) + (x + y)X 2 

+ (1 — x — y — z)(l — x — y)m 2 + z(l — x — y)M 2 . 



(5) 
(0) 
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n 


Numerical Results 


P 


Analytical Results [P=4] 


-30 


-0.150899286980769753D-01 ± 0.771D-26 
0.189229839615898822D-02 ± 0.124D-25 


8 
8 


-0.150899286980482291E-01 
0.189229839615525389E-02 


-80 


-0.405390396284235075D-01 ± 0.580D-15 
0.478581216125478532D-02 ± 0.401D-12 


16 
16 


-0.405390396283445396E-01 
0.478581216112143981E-02 


-120 


-0.608983283726997427D-01 ± 0.556D-15 
0.710062317325786663D-02 ± 0.415D-12 


32 
32 


-0.608983283725815879E-01 
0.710062317309438855E-02 


-160 


-0.812576170752810666D-01 ± 0.549D-10 
0.941543418501223556D-02 ± 0.109D-11 


32 
32 


-0.812576183472699269E-01 
0.941543432496722442E-02 



Table 2: Numerical results of the one-loop vertex diagram with ^fs = 500 GcV, t = — 150 2 
GeV 2 , m = m e = 0.5 x 10~ 3 GcV. Photon mas A = 10™ [GcV] and P denotes P-precision. 
The upper is the result of Real part and the lower is one of Imaginary part. 



s denotes squared central mass energy and m and M are external masses of particles. Again, 
replacing s by s + ie in cq.© and applying extrapolation method, we obtain the numerical 
results. We show the real part of the numerical results with the analytic ones in Table [3] as 
an example [14j . The results are almost consistent each other. However, since high precision 
computation costs huge CPU time, we only have done the quadruple precision arithmetic. 
Thus the results are not reliable when the photon mass is smaller as 10~ 30 GcV. 



4 Dimensional Regularization 



We next consider the diagram contains two massive and two masslcss 
external legs, and internal lines are all massless with the dimensional 
regularization, where D = e + 2 (Fig. [3]). In order to pick up the 
singularities, we apply the sector decomposition technique. The loop 
integral is expressed as follows: 



h = 



(-sxix 3 - 

n=-2,-1.0.— 



tX2X4 — 1113X3X4 — 777,4X1X4) 



2-e 



, (7) 
(8) 



In the case of 



123 GeV 2 ,t = -200 GeV 2 



50 GeV 2 and 



77i| = 60 GeV 2 , the denominator of each integrand corresponding 
to 1/e and the finite term can be zero. Again, +ie is introduced 
in the denominator of eq.([5|) to specify the analyticity of it and the 
extrapolation method does work to perform the numerical integration. 
Table [4] shows the complete numerical integration can reconstruct the 
results with the analytical expressions [17] . 



, ITh 



ITU 




Figure 3: One-loop 
box diagram with 
two massive exter- 
nal legs, so called 
"hard" case. 



5 Summary 

We discussed a new approach of the numerical method for multi-loop integrals. The combi- 
nation of an efficient multi-dimensional integration routine and an extrapolation method is 
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n 


Numerical Results 


P 


Numerical Results [P=4] 


-15 


-0.1927861102670278D-06 ± 0.314D-14 


4,4,2 


-0.1927861122439641E-06 


-20 


-0.2472486348234972D-06 ± 0.586D-15 


4,4,2 


-0.2472486352599175E-06 


-25 


-0.3017111253761463D-06 ± 0.111D-13 


4,4,2 


-0.3017111582758710E-06 


-30 


-0.3562028882831722D-06 ± 0.867D-10 


4,4,2 


-0.3561736812918245E-06 



Table 3: Numerical results of the loop integral of one- loop box diagram with y/s = 500 GeV, 
t = -150 2 GeV 2 , m = 0.5 xl0~ 5 GcV, M = 150 GeV. Photon mas A = 10" [GeV] and P 
denotes P-precision. This is results of Real part. 



n 


Numerical Results 


Analytic Results 


-2 


-0.40650406505E-04 


-0.40650406504E-04 


-1 


-0.34156307031E-03 


-0.34156306995E-03 





-0.14929502492E-02 


-0.14929502456E-02 



Table 4: Numerical v.s. Analytical results in the double precision arithmetic. Real part of 
C n of Fig. [3] with s = 123 GcV 2 , t = -200 GeV 2 , m§ = 50 GeV 2 and mf=60 GeV 2 . 

applicable to the one-loop integral even though it includes infrared singularities both with 
the introduction of the fictitious mass and the dimensional rcgularization. We demonstrated 
the fully numerical calculations for one-loop vertex and box diagram with infrared singu- 
larities. To handle the infrared singularities, high precision control is essential and in some 
case quadruple precision is not enough to obtain the stable results. In our calculations, we 
applied the multi-precision library HMLib which has the high performance of the precision 
control . Our method is efficient and we obtain reliable results for one-loop vertex and box 
diagrams with infrared singularities. 
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